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The one dimensional 5 = 1 XXZ model with next-nearest-neighbor interaction a and Ising- 
type anisotropy A is studied by using a numerical diagonalization technique. We discuss the 
ground state phase diagram of this model numerically by the twisted-boundary-condition 
level spectroscopy method and the phenomenological renormalization group method, and 
analytically by the spin wave theory. We determine the phase boundaries among the XY 
phase, the Haldane phase, the ferromagnetic phase and the Neel phase, and then we confirm 
the universality class. Moreover, we map this model onto the non-linear a model and analyze 
the phase diagram in the a <C — 1 and A ~ 1 region by using the renormalization group 
method. 

KEYWORDS: spin-1, ground state phase diagram, frustration, diagonalization, level spec- 
troscopy, non-linear a model 

1. Introduction 

One dimensional antiferromagnetic Heisenberg spin systems have been an important sub- 
ject in many-body physics for a long time. There are no ordered ground states with a con- 
tinuous symmetry breaking in these systems, because of quantum fluctuation and low dimen- 
sionality. It is well known from Haldane's conjecture 1 ' 2 that there is a fundamental difference 
between half-integer and integer spin chains in their low-lying excitations. In some lattice sys- 
tems, frustration, which is geometrically represented by a triangular lattice, tends to disorder 
the ground state. We expect that competing systems between quantum fluctuation and frus- 
tration plus low dimensionality exhibit fascinating phenomena in their ground states and low- 
lying excitations. Recently, the materials of the S = 1 spin chain with next-nearest-neighbor 
interaction have been reported experimentally (CaV204, 3 NaV(W04)2 4 ). 

In this paper, we study the S = 1 XXZ chain with next-nearest-neighbor (NNN) interac- 
tion, which is described by the following Hamiltonian: 



L 

H = /^(Sf Sf + i + SfSf +1 + A<Sf 5f +1 , 
=i 



+aJ2(SfSf +2 + S? Sf +2 + ASf Sf +2 ) , (1) 
i=i 

where a is next-nearest-neighbor interaction constant and L is system size (L: even). Hereafter, 
periodic boundary condition (PBC) will be assumed unless specifically mentioned. 
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A phase diagram of this model has been studied by the exact diagonalization method for 
the | A| < 2 and < a < 0.3 region, 5 and also by the infinite-system density matrix renormal- 
ization group method (DMRG) for the A > and a > region. 6,7 Many phases of the ground 
state have been found there; the Haldane phase, the double Haldane phase, the Neel phase, 
the double Neel phase, the gapped chiral phase and the gapless chiral phase. The features of 
these phases are explained simply as follows. The Haldane state has an energy gap, where the 
spin correlation decreases exponentially as L — > oo. The Neel state is characterized by twofold 
degenerate ground states, where the staggered magnetization exhibits a long range order. In 
the double Haldane phase, the system can be regarded as two Haldane subchains because 
of competition between the strong next-nearest-neighbor interaction and the weak nearest- 
neighbor interaction. A first-order phase transition occurs between the Haldane phase and 
the double Haldane phase. 8 ' 9 We can also regard the double Neel state as two Neel subchains 
(up-down- up-down or up- up- down- down ) . The gapped chiral state and the gapless chiral state 
have a chiral long range order, whereas the two spin correlation decays exponentially in the 
former state and algebraically in the latter one. 6 ' 7 

On the a = line, it is known that four different ground states exist; the ferromagnetic 
state(A < -1), the XY state(-l < A < 0), the Haldane state(0 < A < 1.17) and the Neel 
state(1.17 > A). 10-13 The ferromagnetic state has a long range order of the magnetization. 
In the XY phase, the spin correlation shows an algebraic decay, and the correlation length is 
infinite. A first-order phase transition occurs between the ferromagnetic phase and the XY 
phase, and the 2D-Ising type phase transition occurs between the Haldane phase and the Neel 
phase. 13 

It has been well known that the Berezinskii-Kosterlitz-Thouless (BKT) transition occurs 
between the XY phase and the Haldane phase, which is explained by the bosonization method 
for a spin-s XXZ model. 14 This transition belongs to the infinite-order phase transition in the 
traditional classification. On this transition point, the logarithmic correction appears in the 
correlation function, the critical exponents and so on. Numerically, the BKT transition has 
been known to occur at the A = and a = point in a spin-1 (S = 1) XXZ model with bond 
alternation constant 5 for \5\ < 0.23, 10,11 and that with uniaxial single-ion-type anisotropy 
constant D for — 2 < D < 0.5. 12 A bosonization study for a spin-s XY model with NNN 
interaction has supported that this particular point is determined as the BKT transition 
point in the 5 = 1 case. 15 Recently, the fact that this point (A = a = 0) is the BKT 
transition point has been analytically proved. 16 For many years, it was a difficult problem 
to determine numerically the critical point and the universality class of the BKT transition, 
because of logarithmic corrections and the slow divergence of the correlation length. Although 
the DMRG method is a powerful method to examine many phases, it is not appropriate 
to determine the BKT transition line. This problem was successfully resolved by the level 
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spectroscopy method, based on the conformal field theory and the renormalization group. 
Moreover, this method was improved by the twisted-boundary-condition level spectroscopy 
method. 18 We use the twisted-boundary-condition level spectroscopy method to determine 
the BKT transition point when A / or a / 0. 

In this work, we present a phase diagram of the S = 1 XXZ chain with NNN interac- 
tion, analyzing the diagonalization data by the twisted-boundary-condition level spectroscopy 
method, 18 the phenomenological renormalization group (PRG) method, 19,20 the spin wave 
theory and the renormalization group method. 

This paper is organized as follows. In the next section, the obtained phase diagram is 
presented. The numerical methods to determine the phase boundaries are explained in §3. In 
§4, we map eq. (1) onto the non-linear a model, and then we analytically examine the phase 
diagram in the region of a <C —1 and A ~ 1. The final section is devoted to a summary and 
a discussion. 

2. Phase diagram 

The obtained phase diagram is summarized in Fig. 1, which consists of the ferromagnetic 
phase, the XY phase, the Haldane phase and the Neel phase. We have determined the BKT 
transition line between the XY and Haldane phases by using the twisted-boundary-condition 
level spectroscopy method. 18 Then, from the analysis of the spin wave theory, we determine 
the first-order phase transition line between the ferromagnetic and XY phases, and that 
between the ferromagnetic and Haldane phases. Moreover, the Haldane-Neel transition line is 
determined by the PRG method. 19 ' 20 

In the a > 0.4 region, it is expected that the gapped and gapless chiral phases, the double 
Haldane phase and the double Neel phase appear in the phase diagram. 6-9 Since the level 
spectroscopy and PRG methods do not work in this region because of incommensurability, we 
do not treat this region in this study. This point will be discussed further in §5. 

3. Numerical methods and results 

3.1 Level spectroscopy and BKT transition 

In this subsection, we briefly review the level spectroscopy method 17 and the twisted- 
boundary-condition level spectroscopy method, 18 and then we determine the BKT transition 
line using the latter method. In addition, we check the consistency of the universality class. 

With the use of the bosonization technique, 14 the one dimensional XXZ spin system can 
be mapped onto the following sine-Gordon model (in Euclidean space-time): 



where <\> = 4>{t, x) is the boson field and a is the lattice constant. After a scaling transformation 
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Fig. 1. Phase diagram of the 5 = 1 XXZ chain with next-nearest-neighbor interaction in the region 
of — 1 < A < 2 and — 1 < a < 0.4. Note that the ferromagnetic stable line for a > 0.4 is given in 
Fig. 6. 



a — > e a, we obtain the following renormalization group equations: 



dK 



-i 



dl 



dj/2 
d^ 



y) y2 - 



(3) 



These equations show that the j/2 term in eq. (2) is an irrelevant field for K > 4 and a relevant 
field for K < 4. On the = 4 line, the BKT transition occurs. Denoting K = 4 + 2yi near 
if = 4, we have 

dj/i(0 /n2 



dl 

d/ 



= -y 2 (/) 2 , 



= -Vi(l)V2(l)- 



(4) 
(5) 



The renormalization flow separatrix is given by y\{l) 2 — 2/2 (0 2 = > 0). These dif- 

ferential equations are the famous Kosterlitz equations. 21 On this separatrix, the solution of 
eqs. (4) and (5) is y2(l) = = 2/1/ (1 + yiZ), where yi is the bare coupling constant and 

I = 1/lnL. 

When the 1/2 term is relevant, the scaling dimension x should be corrected as x — > x + C?/2. 
Since the behavior of some scaling dimensions is different depending on whether in the case 
of relevant 1/2 or in the case of irrelevant 1/2, we expect that the BKT transition point can be 
determined by examining the scaling dimension. 

Now, we introduce the dual field 9 = 9(r, x) to the boson field <\> by the following equations: 



d T (j) = —iKd x 9, d x (j) = iKd T 9. 



(6) 
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The vertex operators are expressed in terms of 9 and <j> as 

V n , m = exp(in\/26>) exp(im-\/20), (7) 

where n and m are integers, conventionally called the electric charge and the magnetic charge, 
respectively. Then, for the free field case (j/2 = 0), the scaling dimensions of the vertex 
operators V njjn are 22 

Xn,m = \{^+ Km2 ) ■ ( 8 ) 

In addition, there is a marginal operator M = (d T 4>) 2 + (d x (j)) 2 , whose scaling dimension is 
given by x marg = 2. From the conformal field theory, the eigenvalues E njTn with PBC are 
related to the scaling dimensions x n ^ m as 23 

E n , m (L) - E g (L) = (9) 

in the limit of L — > oo, where E g (L) is the ground state energy for the size L system under 
PBC, and v s is the spin wave velocity. Using eqs. (8) and (9), we obtain 

^±4,0 = %0,±1 = ^marg (10) 

when K = 4 and 1/2 = (marginal). Considering the invariance under <p <-> —<f>, we obtain 
the combined operators cos(v / 2</ ) ) ~ Vo,+i + ^0,-1 an d sin(\/20) ~ Vb,+i — Vb,-i, whose 
scaling dimensions are xo cos and xosinj 

respectively. In the case of 2/2 / 0, these scaling 
dimensions differ from each other because of the yi cos(\/20) term in eq. (2). Thus, x±4,0j 
^o,cos, ^o,sin and x marg are expressed with the logarithmic correction as x±4,o = 2 — yi(l), 
x ,cos = 2 + 2yi(i)(l + 2i/3), x , si „ = 2 + yi(Z) and x marg = 2 - yi (l)(l + 4t/3), where t 
is a distance from the BKT transition point, defined as 3/2 (0 = + At the BKT 

transition point (t = 0), we have 

X±4,0 = Xmarg- (H) 

Thus we can use this energy level crossing to determine the BKT transition point. 17 

Although the above procedure can be applied to determine the BKT transition point, 
a better method, the twisted-boundary-condition level spectroscopy method, has been pro- 
posed. 18 This method is based on the SU^2)jZi symmetry on the BKT transition point. We 
use the following twisted boundary condition (TBC): 

Sl'Z j = -S*' y ,Si +j = S]. (12) 

Then, the scaling dimensions of the eigenvalues £^~(L) with TBC are given by 



I in 2 / O 2 ' 



or 



x n,m = 2}~K +K \ m+7 >I ' X "' m+1 / 2 ( 13 ) 



E™(L)-E g (L) = 2 -^ X ™. (14) 
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Table I. This table shows relations between quantum numbers (M, P, k) and renormalized scaling 
dimensions x. M — ' ■ BC represents boundary conditions. PBC means periodic boundary 

condition and TBC means twisted boundary condition. Under TBC, the parity is different from 
that of PBC, so we use * to classify. y 2 (0 = ±2/i(0 = + J/iO> where y\ is the bare coupling 

constant and I = InL. 



M 


BC 


P 


k correction in x 


operator in s.G. 


scaling dimension 


±2 




PBC 
TBC 
TBC 


1 

-1* 
1* 


1/2 - yi (/)/4 

1/2 + yi (Z)/4- 3/2(0/2 
l/2 + yi (0/4 + y 2 (0/2 


exp(±i2^/20) 
sin(<^/v / 2) 
cos(<f>/y/2) 


X±2,0 

TBC 
x 0,sin 

TBC 

x 0,cos 



Since we can obtain the half-integer magnetic charge effectively from eq. (13), we now extend 
the magnetic charge m to that including half- integers. From eqs. (8) and (13), we have in the 
case of y 2 = 0, 

X± 2 ,0 = X™°(= X ,±l/2) (15) 

because K = 4. When y 2 / 0, following the above mentioned way, we define the scaling 
dimensions x™ ^ and £™n of the operators cos(0/\/2~) ~ ^o,+i/2 + ^0,-1/2 an d sm((j>/y/2) ~ 
^o,+i/2 ~~ ^b,-i/2) respectively, and also we express these scaling dimensions with the loga- 
rithmic correction as shown in Table I. At the BKT transition point (y 2 = ±yi, yi > 0), we 
have 

^±2,0 = x™c. (16) 

Since the electric charge and the magnetic charge of x± 2 ,o and x™„ are smaller than those of 
x±4,o and x marg , we can obtain more accurate results by using the twisted-boundary-condition 
level spectroscopy method 18 than those by using the original level spectroscopy method. 17 

We can identify the correspondence between the operators in the boson representation and 
the eigenstates of spin chains by comparing their symmetry properties 18 as shown in Table I. 
The Hamiltonian of eq. (1) with PBC is invariant under the spin rotation around the z-axis, 
the translation (S, — > Sj + i) and the space inversion (S, — > Sjv-j+i)- The quantum numbers 
corresponding to these invariance are M = J2i=i^i^ k = 2nn/L(n = 0, • • • ,L/2 — 1) and 
P = ±1, respectively. From Table I, we see that x± 2 ,o corresponds to the lowest excitation 
energy with M = ±2, P = 1, k = 0, and x™n to the lowest excitation energy with M = 0, 
P = — 1* under TBC. The energy level crossing of these two energies for L = 16 and A = —0.2 
is shown in Fig. 2. 

In addition to the logarithmic correction, there is another correction between a lattice 
system and a continuous model. Considering the x = 4 irrelevant field operator L_ 2 Z_ 2 1 
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Fig. 2. Energy level crossing between X2,o and x™^ when L = 16 and A = —0.2. The value of o^ KT 
is obtained to be a? KT = 0.09193. 
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Fig. 3. Extrapolation procedure of af for A = —0.2, where a fitting function a(L) = af + 
Ci/L 2 + C2/L 4 , is assumed. 

in terms of the conformal field theory, 24 we should extrapolate the infinite-size critical point 
o^ KT as 

a(L) = a* KT + C 1 /L 2 + C 2 /L 4 . (17) 

The BKT transition point for A = —0.2 is determined as shown in Fig. 3. The same procedure 
can be carried out varying a with fixed A, and the obtained BKT transition line is represented 
in Fig. 1. 

Next, in order to check the consistency of the universality class, we confirm that two 
universal relations, concerning the central charge c and the energy gap ratio, hold in the 
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Fig. 4. Plot of the effective central charge c as a function of a for A = —0.2. In the massless region, 
c = 1. The BKT transition occurs at a? KT . 



following way. 

The BKT transition line is expected to be described by the conformal field theory with 
the central charge c = 1. It is well known that the finite-size correction to the ground state 
energy is related to the central charge c and the spin wave velocity v s as follows, 23-26 



1 

— E g (L) = £oo - 



7TCV S 



where 



lim J±{E kl (L)-E g (L)}. 



(18) 



(19) 



Here E^ (L) is the energy of the excited state with wave number k\ = ^ and magnetization 
M = 0. In general, we should obtain v s as 

v(L) =v s + d/L 2 + C 2 /L\ (20) 

Although there remain logarithmic corrections in the central charge, which is obtained from 
eq. (18), they are small enough to be 0(l/(ln L) 3 ). Thus we neglect them. It is proper that 
the effective central charge c should correspond to the central charge c in a massless region. 
In a massive region, the correction term in eq. (18) is modified as 1/L 2 — > exp(— L/£). We 
calculate the effective central charge c, the result of which is shown in Fig. 4. We can see that 
the effective central charge changes rapidly from c = 1 (massless) to c = (massive). This 
fact agrees with previous results in refs. 28 and 29. 

Next, we examine the ratio of the energy levels. From eqs. (8) and (9), we obtain 

AE 2 A 



(21) 
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Fig. 5. Plot of AE 2 /AE 1 as a function of a for A = -0.2. In the XY phase, AE 2 /AE 1 = 4, and in 
the gapful phase, AE 2 /AEi = 2. Because of the finite-size effect, AE 2 /AEi varies slowly. 



where AEi = E(M = ±i) - E g (M = 0) (i = 1, 2). This relation is independent of K. In the 
gapful (massive) phase, on the other hand, we have 

(L-oo), (22) 

since the mass of two magnons is twice as large as that of one magnon when each magnons 
are independent. This ratio varying with a in the case of A = —0.2 is shown in Fig. 5. We 
can also see the size dependence of AE2/AE1. Affected by the x = 4 irrelevant field operator 
L-2L-2I, the extrapolated value, described in Fig. 5, can be given as follows, 

= 00 + C 1 /L 2 + C 2 /L 4 , (23) 

where O means AE2 / AE\ . It appears that the ratio is 4 in the massless phase and it changes 
very slowly from 4 to 2 because of the finite-size effect. 
3.2 Instability of ferromagnetic phase 

We determine the unstable boundary of the ferromagnetic phase as follows. If there is no 
intermediate phase, we can obtain the boundary of the ferromagnetic phase by examining the 
energy level crossing between the energy of the ferromagnetic state, E = (1 + a)AL, and that 
of the M = state. Moreover, we compare the numerical result with that obtained by the 
ferromagnetic spin wave theory (see APPENDIX A). 

The phase diagram near the ferromagnetic phase determined by these two methods is 
shown in Fig. 6. Below a = 0.38, the ferromagnetic stable line obtained from the diago- 
nalization data agrees with the line determined by the spin wave theory. The fact that the 
convergence of the ferromagnetic boundary depending on L is not so good above a = 0.38 
suggests that the incommensurability begins to modulate the ground state energy. The be- 
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Fig. 6. Phase boundary of the ferromagnetic phase, obtained by the energy level crossing from the 
diagonalization data and the spin wave theory. 



havior of the ferromagnetic phase boundary obtained by the energy level crossing is, however, 
analogous to that by the spin wave theory. Therefore, we adopt the line of the spin wave 
theory as the ferromagnetic phase boundary in Fig. 1. 
3.3 Haldane-Neel transition 

The phase transition between the Haldane phase and the Neel phase belongs to the second 
order phase transition. In order to determine the transition line between these phases, we 
apply the PRG method. 19 ' 20 This method is explained as follows. In the Neel phase, the 
ground state is twofold degenerate in the thermodynamic limit. For finite L, however, this 
degeneracy is splitting into two energy levels. The energy gap between these energy levels is 
AE(L,a,A) = E(M = 0, P = -1, k = tt) - E g (M = 0, P = 1, k = 0), which behaves as 
AE ~ exp(— L/£) in the L — > oo limit. On the other hand, the Haldane state has the energy 
gap, which remains finite in the thermodynamic limit. Thus the product LAE(L) increases 
(decreases) with L for A < A™ (A > A^ N ), where A^ N means the critical value of A 
at the Haldane-Neel transition. Therefore, we use the following equation to determine the 
Haldane-Neel transition: 

LAE(L, a, A C (L, L + 2)) = (L + 2)AE(L + 2, a, A C (L, L + 2)). (24) 

The A-dependences of LAE(L) with a = —0.5 for various L are shown in Fig. 7. We extrap- 
olate the infinite-size critical point A^ N as follows: 

Ac(i , i + 2HAc HN + _£_ + _^_ (25) 

When a = —0.5, the infinite-size transition point A^ N is determined as shown in Fig. 8. The 
same procedure can be carried out by varying A with fixed a. The obtained Haldane-Neel 
transition line is shown in Fig. 1 and also in Fig. 9 together with the finite-size results. 
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Fig. 7. The A dependence of LAE(L) with a = -0.5 for L =10, 12, 14, 16 and 18. The crossing 
points (•) are the finite-size transition points. 
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From the quantum Monte Carlo simulation in the case of a = 0, 13 this transition has 
been known to belong to the 2D-Ising type universality class. Therefore, we expect that the 
universality class of the whole Haldane-Neel transition line is of the 2D-Ising type. 

4. The analysis of the anisotropic non-linear er model 

In order to examine whether the Haldane phase does or does not survive in the region of 
a <C — 1 and A ~ 1, we perform the following analytical calculation, the details of which is 
discussed in APPENDIX B. First, we map the spin-s (s is integer) XXZ spin chain with NNN 
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Fig. 9. Phase diagram near the Haldane-Neel transition line. 



interaction onto the non-linear a model with Ising-type anisotropy A in the large-s limit. 
Then, we obtain the renormalization group equation (eq. (B-19)) from the non-linear a model 
by using the renormalization group method. 33 From this equation, the renormalization group 
flows are represented by A and a as 

| A - 1| oc ^ 1 ~ 4a exp (-271-V1 - 4a), when |A-1|<1. (26) 
1 — a 

We assume that these renormalization group flows are connected to each renormalization 
group fixed point: the XY or the Ising fixed point. These fixed points may be asymmetric on 
the phase diagram. Therefore, we can draw the schematic phase diagram, comparing with the 
numerical data, as shown in Fig. 10. The numerical results are in good agreement with the 
analytical lines near A ~ 1 and a ~ — 1. 

The phase is massive at the isotropic point, A = 1 and a — > — oo, since the mass term (eq. 
(B-20)) exists. Thus, we may conclude that the phase transition between the XY phase and 
the Neel phase does not directly occur in this limit and the Haldane phase appears between 
these two phases. 

5. SUMMARY AND DISCUSSION 

The ground state phase diagram of the S = 1 XXZ chain with NNN interaction has been 
determined numerically by analyzing the numerical diagonalization data using the twisted- 
boundary-condition level spectroscopy method 18 and the phenomenological renormalization 
group method, 19 ' 20 and analytically by the spin wave theory. Particularly, it has been re- 
solved that the phase transition between the XY phase and the Neel phase does not occur 
directly for a <C —1 near A = 1 by mapping onto the non- linear a model 30-32 and using the 
renormalization group method. 33 

In this paper, we have not treated the region for a > 0.4, except for the vicinity of the 
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Fig. 10. Phase diagram near the a = — 1 and A = 1 point. We compare the schematic phase diagram 
obtained from eq. (26) with the numerically extrapolated data to L — > oo. The numerical results 
are in good agreement with the analytical line near A = 1 and a = — 1. 

ferromagnetic phase, because of the incommensurability. When A = 1.0, the commensurate- 
incommensurate (C-IC) change is known to occur at a = ac-ic, where ac-ic ~ 0.28 is 
obtained by the DMRG method 8,9 and ccc-ic ~ 0.38 by the diagonalization technique. 34 
Note that this C-IC change is not a phase transition because the energy gap does not vanish. 

In principle, the phase transition line between the Haldane phase and the double Hal- 
dane phase can be determined by Kitazawa's twisted-boundary-condition level spectroscopy 
method. 35 However, we failed to determine this transition line by this method. Since the 
ground state energy is modulated with the long wave length in the incommensurate region, 
the extrapolation procedure is difficult. In this region, some long range orders, chiral order 
and string order, exist, and therefore the infinite-system DMRG method is known to be rel- 
atively useful. 6 ' 7 However, it is difficult to determine phase boundaries of the gapless chiral 
phase, where the string order parameter is exhibiting a quasi-long-range order, by using the 
DMRG method. In order to advance the analysis of the incommensurate region, we need to 
investigate the incommensurability and we must find an effective method to treat this nasty 
problem. 
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Appendix A: The analysis of spin wave instability 

In this Appendix, we consider the instability analysis of the ferromagnetic phase of the 
spin-s XXZ chain with NNN interaction. We consider the Hamiltonian of eq. (1), but S x ,S y 
and S z are spin-s operators now. When we denote the fully ferromagnetic state by |0), we 
obtain the ground state energy for the fully ferromagnetic phase as 

n\0)=E \0), (Al) 

where Eq = As 2 L + aAs 2 L and s is the eigen value of S z for |0). Then, in order to see the 
lowest excitation, we consider the following operator which is the linear combination of the 
spin lowering operator, 

J2S m exp(-ikm), (A-2) 

m 

where k is the wave number. When this operator acts on |0), we obtain 

ftJ]S-exp(-ifcm)|0> = {E + uj{k)}J2s-exp{-ikm)\0), (A-3) 

m m 

where u(k) represents the dispersion relation, 

uj(k) = 2s{(cos k — A) + a(cos 2k - A)} (A-4) 

= ^{ 2 «( cosfc + i) 2 -^-«- A -« A }- ( A - 5 ) 

If oj(k) is negative, the fully ferromagnetic state is unstable because of the spin wave excitation. 
Therefore we can regard the line of uo{k) = as the ferromagnetic stable-unstable boundary. 
Since the prefactor 2s is positive, we omit 2s in the following analysis. 

From eq. (A-5), we can see that the minimum value of the dispersion \uj(k)\ is different 
whether in the case of a > 1/4 or in the case of a < 1/4. In order to examine the incommen- 
surate region, it is enough to consider for a > 0. 

In the former case, \uj(k)\ has a minimum value when k — i arccos Then, we have 

A = - 1 + 8a \- (A-6) 
8a (1 + a) V ; 

Therefore the fully ferromagnetic phase becomes unstable for k = ± arccos ^ modes of the 

spin wave excitation. Classically, these modes are connected with incommensurability. 

In the latter case, |oj (fc) | is minimum when k = ±ir. Then, we have 

a — 1 , . 

A = -— - (A-7) 
a + 1 

Therefore the fully ferromagnetic phase becomes unstable for k = ±tt mode of the spin wave 
excitation. 
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Appendix B: The non-linear <x model 

The spin-s XXZ chain can be mapped onto the non-linear a model with Ising-type 
anisotropy in the large-s limit. 30-32 

We first consider the antiferromagnetic Heisenberg Hamiltonian with NNN interaction a: 

H = J2 S i + (B-l) 

i i 

with a < 0. We want to keep low energy modes as follows, 

Si-i-lYspi + ali, (B-2) 

where cp corresponds to the Fourier mode with the wave number k ~ ir, and 1 corresponds to 
the Fourier mode with the wave number k ~ 0. These operators are defined as 

1\ (S2i+1 — S2i) 



if 2i 



l[2i+ 1 - 



2s 

(Sgi+l + S2i) 

2a ' 



(B-3) 



(B-4) 



where a is the lattice spacing. Note that ip and 1 are slowly varying function of i. In the 
continuous limit, these operators obey the following commutators, 

>W(y)] = ie a ^P(x)S(x-y), (B5) 

r(x),/(y)] = ie a ^^(x)5(x-y), (B-6) 

99 a (x),/(y)] = ie a ^r(x)5(x - y)2a^ ^ 0, (B-7) 

where {a, (3, 7} = {x, y, z) and e a/?7 is the three-dimensional Levi-Civita symbol. Repeated 
indices are normally summed. These operators satisfy the following constraints, 

<p-\ = 0, (B-8) 

M 2 = 1- (B-9) 
Now, the Hamiltonian density can be mapped onto the 0(3) non-linear a model, 



' 1 ' An dx) g \dx 



(B-10) 



where the velocity i> = 2sa, the topological angle 6* = 2-7rs and the coupling constant g = 
2 

Next, we consider the Hamiltonian with Ising-type anisotropy A as follows, 

H = • S m + a S t ■ S t+2 + £ {(A - 1) SfSf +l + a (A - 1) SfSf +2 } . (B-ll) 



We only choose the leading term in the anisotropic terms as 



SfS. 



S; S, 



i u i+2 



2 2 



2 2 

S <p z . 



(B12) 
(B-13) 
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Then, we obtain the following effective Hamiltonian density, 

6 d(p\ 2 1 / cfy?\ 2 l vm 2 2 




"^i'i'-tij^iijrv" (B44) 

where 

m2 _ 2(A-D(-l + o) 
a 2 (l -4a) 1/2 

We have obtained the Hamiltonian density of the anisotropic non-linear a model. The corre- 
sponding Euclidean (2-D classical) Lagrangian density is 30-32 

C = Yg { {8 ^ )2 + ^ + L^* ' X dv(p) ] (B ' 16) 

where we set v = 1 for this discussion. Since we take an integer spin case, we can neglect the 
topological term. 

Renormalization group equations of this model without the topological term are obtained 
by Nelson and Pelcovits. 33 Since d = 2 (dimension) and n = 3 (number of components) in 
our case, we obtain the differential equations for the "dressed" coupling constant g and the 
"dressed" mass rh: 

d ~ 9 - 1 ( ~ g2 ^ (B-17) 



dln& 2vr \\ + rh 2 

*^ = (B.18) 
d In b 7r 1 + m z 

Removing In 6, we can obtain the renormalization relation between rh 2 and g, 

dm 2 47r o /-, ~2\ 2 "t- 2 /t^ ^\ 

_ = m 2 h +m 2 ) . (B-19) 

dg g 2 > g 

Therefore, the relation between the dressed rh 2 and the dressed g is 

m 2 (g) oc exp ( J . (B-20) 



9 J 

For rh 2 > 0, all flow lines terminate in a fixed line at rh 2 = +oo. One of these lines is believed 
to terminate at the BKT transition point g c . Therefore, this flow line should be a critical line 
between the XY phase and the disordered (Haldane) phase. From the similar discussion, we 
obtain the Ising critical line for rh 2 < case. Finally, in terms of the original Hamiltonian 
(eq. (B-ll)) by the use of eq. (B-15), we now obtain the critical line eq.(26) expressed as 



I A - II oc — — — exp (-2W1 - 4a). (B-21) 
I — a 

Note that eq. (B-21), or eq. (26), can be extended to the a — > — oo limit. 
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